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Abstract 

Four scenarios have been proposed for the low-temperature phase behavior of hquid water, each 
predicting different thermodynamics. The physical mechanism which leads to each is debated. 
Moreover, it is still unclear which of the scenarios best describes water, as there is no definitive 
experimental test. Here we address both open issues within the framework of a microscopic cell 
model by performing a study combining mean field calculations and Monte Carlo simulations. We 
show that a common physical mechanism underlies each of the four scenarios, and that two key 
physical quantities determine which of the four scenarios describes water: (i) the strength of the 
directional component of the hydrogen bond and (ii) the strength of the cooperative component of 
the hydrogen bond. The four scenarios may be mapped in the space of these two quantities. We 
argue that our conclusions are model- independent. Using estimates from experimental data for H 
bond properties the model predicts that the low-temperature phase diagram of water exhibits a 
liquid-liquid critical point at positive pressure. 

PACS numbers: 
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Water's phase diagram is rich and complex: more than sixteen crystaUine phases and 

tinn 

two or more glasses [2, 13|, 1^ have been reported. The liquid state also displays interesting 
behavior, such as the density maximum for 1 atm at 4°C. The volume fluctuations {{SV)"^), 
entropy fluctuations {{6S)'^), and cross-fluctuations between volume and entropy {6V6S), 
proportional to the magnitude of isothermal compressibility Kt, isobaric speciflc heat Cp, 
and isobaric thermal expansivity ap, respectively, show anomalous increases in magnitude 
upo„ cooliug B. Further, these quantities display au appareut d.vergeuce for 1 at,„ at 
— 45°C [5[, hinting at interesting phase behavior in the supercooled region. 

Microscopically, water's anomalous liquid behavior is understood as resulting from the 
tendency of neighboring molecules to form hydrogen (H) bonds upon cooling, with a decrease 
of local potential energy, decrease of local entropy, and increase of local volume due to the 
formation of local open structures of bonded molecules. Different models include these 
H-bond features, but depending on the assumptions and approximations of each model, 
different conclusions are obtained for the low-T phase behavior. The relevant region of the 
bulk liquid state cannot be probed experimentally, and none of the theories tested, because 
crystallization of bulk water is unavoidable below the homogeneous nucleation temperature 
Th (-38°C at 1 atm). 



FOUR SCENARIOS FOR SUPERCOOLED WATER 



Due to the difficulty of obtaining experimental evidence, theoretical and numerical analy- 
ses are useful. Four separate scenarios for the pressure-temperature (P — T) phase diagram 
have been proposed: ^ 

(1) The stability limit (SL) scenario [6] hypothesizes that the superheated liquid-gas spin- 
odal at negative pressure re-enters the positive P region below Th{P)- In this view, the 
liquid state is delimited by a single continuous locus, Ps(T), bounding the superheated, 
stretched and supercooled states. There is no reference to the phase into which the liquid 
transforms when P Ps{T). As the spinodal is approached, Kt, Cp, and ap oo. A 
thermodynamic consequence of the SL scenario is that the intersection of the retracing spin- 
odal with the liquid-vapor coexistence line must be a critical point [2]. The presence of a 
such a critical point in the liquid-vapor transition, although possible, is not confirmed by 
any experiment. This fact poses a serious challenge to the SL scenario. 
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(II) The liquid-liquid critical point (LLCP) scenario [3] hypothesizes a first-order phase 
transition line between two liquids — a low density liquid (LDL), and a high density liquid 
(HDL) — which terminates at a liquid-liquid critical point C. HDL is a dense liquid with a 
highly disordered structure, whereas LDL has a lower density and locally tetrahedral order. 
The experimentally observable high density amorphous (HDA) and low density amorphous 
(LDA) solids correspond, in this scenario, to a structurally arrested state of HDL and LDL 
respectively . Starting from C, the locus of maxima of the correlation length C, (the 
Widom line) projects into the one-phase region 10|. Asymptotically close to the critical 
point, response functions can be expressed in terms of ^, hence, these too will show maxima, 
e.g., as a function of T upon isobaric cooling. These maxima will diverge upon approaching 
C Furthermore, for P > Pc, the pressure of C, the response functions will diverge by 
approaching the spinodal converging to C. Specific models suggest [7|, that Pc > 0, but 
the possibility Pc < has also been proposed 12l |. 

(HI) The singularity-free (SF) scenario [l3j hypothesizes that the low-T anticorrela- 
tion between volume and entropy is sufficient to cause the response functions to increase 
upon cooling and display maxima at non-zero T, without reference to any singular behav- 
ior. Specifically, Sastry et al. [13] consider the temperature of maximum density (TMD) 
line, where the density has a maximum as a function of temperature, and prove a general 
thermodynamic theorem establishing the proportionality between the slope of the TMD, 
(dP/dT) 

TMD, and the temperature derivative of Kt- Thus, since the TMD has negative 
slope in water, it follows that (c?P/c}T)tmd < 0, and therefore Kt must increase upon 
cooling, whether there exists a singularity or not. 

(IV) The critical-point free (CPF) scenario [l^ hypothesizes an order-disorder transi- 
tion, with possibly a weak first-order transition character, separating two liquid phases and 
extending to P < down to the superheated limit of stability of liquid water. This scenario 
effectively predicts a continuous locus of stability limit spanning the superheated, stretched 
and supercooled state, because the spinodal associated with the first-order transition will 
intersect the liquid-gas spinodal at negative pressure. No critical point is present in this 
scenario. 

These four scenarios predict fundamentally different behavior, though each has been 
rationalized as a consequence of the same microscopic interaction: the H bond. A question 
that naturally arises is whether the macroscopic thermodynamic descriptions are in fact 
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connected in some way. Previous works have attempted to uncover relations between several 

3 or (II) and (III) 16, 17|. Here 
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of the scenarios, for example between (I) and (II) 
we offer a relation linking all four scenarios showing that (a) all four can be included in 
one general scheme, (b) the balance between the energies of two components of the H bond 
interaction determines which scenario is valid. Morevover, we argue that current values for 
these energies support the LLCP scenario. 



COOPERATIVE CELL MODEL OF WATER 



We analyze a microscopic model [18| of water in which the fluid is divided into cells 
with nearest neighbor (n.n.) interactions. The division is such that each cell is in contact 
with four n.n., mimicking the first shell of liquid water. The case of a shared H bond, due to 
more than four molecules in the first shell, is assimilated with the case in which a H bond is 
broken, since the interaction energy of a shared bond is less than half the energy of a single 
H bond Qiao]. 

The goal of the model is to represent, microscopically, the essential features of the inter- 
action among water molecules, while being able to qualitatively understand the importance 
of each of these features. To this end the interaction among cells is separated into four 
distinct components. 

The first component of the interaction is due to the short-range repulsion of the electron 
clouds. This is incorporated into the model by assigning to each cell i G [1, A^] (a) a volume 
Vi > fo, where Vq is the exclusion volume per molecule, and (b) a maximum of one molecule. 

The second component includes all the isotropic long-range attractive interactions, such 
as the instantaneous induced dipole-dipole (London) interactions between the electron clouds 
of different molecules or the isotropic part of the hydrogen bond [2?]. We refer to this 
component as the van der Waals attractive interaction, keeping in mind, however, that 
this component includes not only the (weak) London dispersion interaction, but also the 
(stronger) isotropic interaction of the hydrogen bond. The overall sum of the isotropic — 
attractive and repulsive — interactions can be represented in different ways. The one we 
adopt in a mean field (MF) treatment is 

J^ = -e^ninj, (1) 
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where if v^/vi < 0.5 we set the index rii = 0, and if 0.5 < Vo/vi < 1 we set Ui = 1, hence 
nj = if the density in the cell is gas-like, and rij = 1 if the density in the cell is liquid-like; 
e > is the characteristic energy of the attraction and the sum is over all n.n. pairs {i,j)- 

The characteristic feature of H2O is its ability to form H bonds between neighboring 
molecules. This interaction has a strong directional component due to the dipole-dipole 
interaction between the highly concentrated positive charge on each H and each of the two 
excess negative charges concentrated on the O of another water molecule. Accordingly, 
the third component incorporated here is this orientational-dependent interaction, which 
includes the covalent component of the bond [22]. To account for the orientational degrees 
of freedom of each water molecule, we assign to each cell i four bond variables aij = 1, ...,q 
(one for each n.n. cell j), representing the orientation of molecule i with respect to molecule 
j. We choose the parameter q by selecting 30° as the maximum deviation from a linear 
bond, i.e. q = 180°/30° = 6, hence every molecule has = 6^ = 1296 possible orientations. 



(The effect of choosing a different value for q has been analyzed in 23||.) We say that a bond 
is formed between cells i and j if a^j = aji. 

Experiments show that formation of the H bonds leads to an open — locally tetrahedral — 



2J]. This effect is incorporated 



structure that induces an increase of volume per molecule 
in the model by considering the total volume to be given as 

V = Nvo + NubVub, (2) 

where 

NiiB = ^ninj5^^^^^^. (3) 

is the total number of H bonds, 6a^b = 1 if a = Sa^b = if a 7^ 6, and f hb is the volume 
increase per H bond [isl. Bond formation also leads to a decrease in the local potential 
energy, hence we add to the Hamiltonian in Eq. ([1]) the term 

^HB = -JiVHB, (4) 

where J > is the characteristic energy of the directional component of the H bond. 

Another key experimental fact is that at low T the 0-0-0 angle distribution in water 



becomes sharper around the tetrahedral value 251], suggesting an interaction that induces 



a cooperative behavior among bonds. For water, four-body and higher order interactions 
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seem to be negligible with respect to the three-body term [26|, |27|. Hence, the fourth 
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component to the interaction potential is the many-body effect due to H bonds 
which minimizes the energy when the H bonds of nearby molecules assume a tetrahedral 
orientation. This is accomplished by further adding to the Hamiltonian in Eqs. ([T]) and 
the term 

'^oop = ~Ja ^ ] ^ ] '^o-ifc,(Ti£, (5) 

where > is the characteristic energy of the cooperative component of the H bond, and 
{k,i)i indicates one of the six different pairs of the four bond variables of molecule i. This 
interaction introduces a cooperative behavior among bonds, which may be fine tuned by 
changing J„. Choosing Jcr = leads to H bonds which form independent of neighboring 



bonds 



13|, while J(j ^ oo leads to fully dependent bonds 



now given by 



coop • 



3l| . The total Hamiltonian is 

(6) 
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his model is studied using both MF analysis and Monte Carlo (MC) simulations 23|, 



33, 



34, 



Details of the MF and MC techniques are available elsewhere [23 



36j. In 



the following we adopt J = J/e, J„ = Ja/t and ^hb = 2fo. 



MEAN-FIELD RESULTS 



Three qualitatively different phase diagrams are found, dependent on the strengths of the 
H bond energy parameters, J and J„ (Fig. [1]). 

When = the model coincides with that proposed in [l^, which gives rise to the SF 
scenario (Fig. [T^). For Q < J„ < J 12 the model displays a liquid-liquid transition ending 
in a LLCP at Pc > (Fig. [Do) [3^. For J/2 < J^/ < a + bJ, where a = 0.30 ± 001 and 
b = 0.36 ±0.01 are fitting parameters, a LLCP occurs at Pc < (Fig. [It). For > a + bJ, 
a liquid-liquid transition with no critical point is found, consistent with the CPF scenario 
(Fig. [T]i). In Fig. [2] we summarize these results in the J vs. parameter space. 



Limiting behavior between the four cases 



In the following we discuss how, by tuning J and J^, we can pass from one scenario to 
another in a continuous way. 



(i) By beginning with the LLCP scenario, and studying the hmit — * 0, we find Tc 0. 
Moreover, we find that Kjn and ap diverge as |T — Tc'\^^ for any value of J„, including 
Jo- ^ and Tc 0. Further, we find for the entropy S that, for any value of J^, 
{dS/dT)p oc |T - Tc'\-^. Hence Cp = T{dS/dT)p diverges as |T - Tct^ when Tc > 0. 
When Tc' = {J„ = 0), Cp is constant, as in the SF scenario [37|. Therefore, the SF 
scenario coincides with the LLCP scenario in the limiting case of Tq' 0, for J^- ^ (Fig. 

(ii) Again, beginning with the LLCP scenario, and increasing J„ while keeping other 
parameters constant, we observe that C moves to larger T and lower P, with Pc < for 

Jo > J/2 (FigHb). 

(iii) With further increase of J„, C approaches, and eventually reaches, the liquid-gas 
spinodal. For larger values of J^ only the liquid-liquid transition remains, which is precisely 
the CPF scenario [l4| (Fig. [T}i). Hence the CPF scenario differs from the LLCP scenario 
only in that C is now inaccessible, lying beyond the region of liquid states. The same result 
may be obtained by decreasing J, while fixing J^ and other parameters. Here a decrease of 
J moves C to lower pressure, i.e. towards the liquid spinodal, while the entire liquid-liquid 
phenomena moves to succesively lower temperature. In all cases, the location of C varies 
continuously with variation of J and J^. 

(iv) In the case of the CPF scenario, we find that the su per heated liquid-gas spinodal 



merges with the supercooled liquid-liquid spinodal, as in Ref. Q. This gives rise to a liquid 
spinodal which retraces in the P-T plane. This feature resembles the main characteristic 
of the SL scenario, where the high-T liquid has a limit of stability at P < that retraces 
toward P > at low T. Here this retracing locus is formed by two spinodal lines, with 
different signs of the slope, that merge at P < 0. Therefore, in the framework of the present 
model, the CPF scenario and the SL scenario coincide, corresponding to the case in which 
the cooperative behavior is very strong. 

Linearity of the lines separating one scenario from another in J Jo plane 

For the cell model, we can derive 

Tc = Jja + ff{Jl) (7) 
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and 

Pc = (r/vuB) + PTc + ^(T^,). (8) 

Here a > and (3 < are constants and, in the MF context, J* = J + 3J„. Symbols 
where X is J„ or Tq', represent terms of order or higher, that are neghgible 
when X <^ 1. Our MF results confirm the relations in Eq. ([7j) and ([H]), with a ~ 0.74/cs/e 
and (3 ~ —TAks/vo, with negligible (^(X^) terms. 

Therefore, we can rewrite the above relations as J — Pc'Vub/^ = — (3 + /9t'HB/tt) Jo- = 2Ja-, 
when fHB = 2fo. As a consequence, for the case Pc = 0, we find = J/2, which is exactly 
what we find numerically in Fig. [2] along the line separating the LLCP scenario with Pc > 
(valid for J„ < J/ 2) and the LLCP scenario with Pc' < (valid for J„ > J 12). 

It is possible to show that Eq. ^ can be generalized to Pll = {J* /vub) + PTll + ^(Tlj^), 
where Tll and Pll are the T and P along the liquid-liquid transition line. Our MF results 
are in good agreement with this prediction. 

We can estimate the equation of the line separating the LLCP scenario with Pc < and 
the CPF/SL scenario in the J- Jo- plane, by using the Eq. (IHl), together with the equation for 
the liquid-gas spinodal. In particular, we adopt a parametric fit, in terms of the parameter 
J, of the spinodal pressure with respect to the spinodal temperature, and we evaluate the 
line separating the LLCP and CPF/SL scenarios for J ^ when C is on the spinodal. 
From this approximate approach, we derive that = + 7 J, with ~ 0.2, of the same 
order of magnitude of the fitting parameter a ^ 0.30 in Fig. [2l Yet, 7 7^ 0.36, the value of h 
in Fig. [2], as a consequence of the strong approximations made. 



MONTE CARLO RESULTS 



To test the validity of our MF calculations, we perform MC simulations in the NPT 
ensemble 361]. To this end, 

(i) we consider that the total volume is = Vmc + -^HB'i'HB, where Vmc ^ -^"^o is a 
dynamical continuous variable; 

(ii) we assume that the system is homogeneous with all the variables rii set to 1; with 
this assumption the gas state occurs when p = N/V < 0.5/fo; 

(iii) we replace the isotropic repulsive and attractive terms of the Hamiltonian in Eq. ([6]) 
with a Lennard- Jones potential, more suitable for continuous distances r between particles, 
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with attractive energy e > plus a hard-core repulsion at distance tq 



oo if r ^ ro, 

Uw{r)^{ , (9) 



if r > vq. 



Here tq = (vq)^^^ and ci is the system dimension [33] (the hard-core repulsion reduces the 
computational cost and does not change the phase diagram); the distance between two n.n. 
molecules is (V/NY^'^, and the distance r between two generic molecules is the Cartesian 
distance between the centers of the cells in which they are enclosed. 

(iv) We consider the system in ci = 2 dimensions. While the MF results are valid for any 
dimension so long as the number of n.n. molecules is four, the MC results hold for a system 
with coordination number four and two dimensions. Since the results in the two cases are 
qualitatively comparable, we do not expect a strong dependence of the phase diagram on 
dimension. 

We simulate this system for = 10'' molecules arranged on a square lattice, adopting 
Wolff's algorithm to equilibrate at low T [36], for different values of J^, keeping constant 
J = 0.5, and vub/vq = 0.5 (Fig. E]). 

For large values of J„ {J„ = 0.5 > a + bJ), we find a HDL-LDL first-order phase 
transition that merges with the superheated liquid spinodal as in the CPF scenario (Fig. [3^). 
At lower J„ {J„ = 0.3 > J/2), a HDL-LDL critical point appears at P < 0, from which 
emanates the locus of Cp maxima (used here as an approximation of the liquid-liquid Widom 
line), which intersects the superheated liquid spinodal (Fig. [3b). By further decreasing J^r 
{Jfj = 0.05 < t//2), the HDL-LDL critical point occurs at P > 0, with the line of Cp 
maxima intersecting the P = axis (Fig. [St). For J„ = 0.02, approaching zero, we find that 
the temperature of the HDL-LDL critical point approaches zero and the critical pressure 
increases toward the value P = e/vo independent of Ja- In this case, we can show that 
Eq. (IH]) still holds, but with J* = J. The line of Cp maxima approaches the T = axis for 
Jcr ^ 0. These results confirm the qualitative behavior found with the MF calculations. 



COMPARISON WITH OTHER THERMODYNAMIC MODELS 

To show that our analysis offers a general framework within which to analyze the super- 
cooled water phase diagram in terms of the interplay between the strengths of the directional 
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contribution to the H bond interaction and its cooperative part, we compare our results with 
those from other thermodynamic models that can reproduce more than one scenario by tun- 
ing appropriate parameters jll, 



15 



17|. 



One free energy model with cooperative interactions is the one introduced by Tanaka 



17(1 . He shows that, as in the SF scenario, water's anomalies are the effect of the excitation 



of locally favored structures upon cooling, which have lower energy and larger volume than 
normal-liquid structures. As in our model, in Tanaka's model increasing the cooperativity 
among excitations of locally favored structures leads to the LLCP scenario. Moreover, 
Tanaka's model LLCP is regulated by relations such as our Eq.(I71) and (IHl). Therefore, by 
increasing the strength of the cooperative interaction, the LLCP will eventually reach the 
limit of stability of the liquid, as in the CPF/SL scenario. 



We next consider the free energy model introduced by Poole et al. [ll|, in which a van 
der Waals free energy is augmented to include the effect of H bond formation. The H bond 
interaction is characterized by two free parameters: the strength of the H bond, and a 
geometrical constraint on H bond formation. The fraction of molecules that form H bonds 
with decreased energy and entropy is determined by a distribution over molar volumes, the 
width of which is A. Poole et al. show that, by keeping A fixed, their model displays a SL 
scenario for weak H bond energy, and a LLCP at positive pressure for strong H bond energy. 
This corresponds in our model to increase the H bond coupling J from J < {J„ — a)/b to 
J > 2 Jo-, while keeping J^j > a fixed. 

Next we study the effect of varying the other H bond parameter in the Poole et al. model, 
the width A. Keeping the H bond energy fixed, we produce the LLCP phase behavior at 
large A and the SL phase behavior at small A. Hence a decrease of A has the same effect 
on the phase diagram as an increase in the H bond cooperativity in our model. 

This result is consistent with that of Borick et al. [15] for their Hamiltonian model that 
incorporates the cooperativity of H bonds trough the same mechanism used by Poole et 
al., i.e. by adopting a distribution with width A that makes the H bond strength density 
dependent. By decreasing A, Borick et al. find that the LLCP moves to lower P and higher 
T. This behavior makes sense physically, as a more all-or-nothing distribution of H bonds 
(small A) implies a more cooperative process of bond formation. It also implies that the 
models of Poole et al. and Borick et al. give rise to the SF scenario only in the limiting case 
of infinite A. 
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We conclude that all four models give a consistent physical picture. This suggests that 
our result, expressed in terms of strength of the directional and cooperative components of 
the H bond, as summarized in Fig. [21 is general. 



ESTIMATES FROM EXPERIMENTAL DATA 



In the framework of the scheme presented here, in which directionality and cooperativity 
are the two relvant physical parameters, we propose that the way to understand which 
scenario best describes water is to probe the energy of the covalent part of the H bond 
interaction 



22l | and the energy of the cooperative component of the H bond interaction 



281, 



29l . |30|. Experiments measure H bonds in ice Ih to be approximately 3 kJ/mol stronger 



than in liquid water 



381]. Attributing this increase to a cooperative interaction among H 



bonds [39|, we can estimate the value of in the cell model to be ~ 1.0 kJ/mol. An estimate 
of the van der Waals attraction, based on isoelectronic molecules at optimal separation, 
yields e ~ 5.5 kJ/mol 40|]. The optimal H bond energy, -Ehb, has been measured to be 
:^23.3 kJ/mol 4l(]. By considering tetrahedral clusters of H bonded molecules, with H bond 
and van der Waals interactions between n.n. molecules (and appropriately reduced van 
der Waals interactions between second and third n.n. molecules), we derive the value for 
the directional component of the H bond, J ~ 12.0 kJ/mol. Other experimental estimates 
suggest that breaking the directional component of the H bond requires ~ 6.3 kJ/mol 42 1. 

Both estimates from experiments fall within the range of 1.1 < J < 2.2, with J^- ~ 0.2, 
i. e. with Jcr < ^/2. Therefore, within our model, these values lead to the LLCP scenario 
with Pc > 0. In particular, MF calculations with = 0.2, J = 2.2 and whb = 2fo, predict 
a LLCP at Tc = 0.25e/kB and Pc = ?>.he/ {v^kB). 



CONCLUSIONS 



We have shown that a microscopic cell model of water, by taking into account the co- 
operativity among H bonds, is able to produce phase behaviors consistent with any of the 
proposed scenarios for water's phase diagram. It is the amount of cooperativity in rela- 
tion to the strength of the directional component of the H bond that establishes which 
scenario holds. For no amount of cooperativity, the SF scenario is recovered. By increasing 
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the amount of cooperativity in relation to the H bond directional strength, a liquid-liquid 
transition grows out from the T = axis, ending in a LLCP. With sufficiently strong cooper- 
ativity, this LLCP lies beyond the region of stable liquid states, leaving only the liquid-liquid 
transition, consistent with the CPF scenario. In this case the spinodal associated with the 
transition acts as the line predicted in the SL scenario. 

Comparison with previous models gives consistent results. Hence we argue that each 
of the four scenarios proposed for the phase diagram of liquid water may be viewed as a 
special case of our general scheme. This scheme is based on the assumption that water- 
water interaction is characterized by an isotropic component, a directional component and 
a cooperative component, and that H bond formation leads to an open local structure. 



43 



44 



45 



, |46|, |47[ or only on 



Alternative mechanisms, based only on isotropic interactions 
directional interactions 



48] have been considered and their relevance for the water case is 
an open question. Finally, estimates for the three components of the H bond interaction, 
based on experimental data, lead to the conclusion that the LLCP scenario with a positive 
critical pressure holds for water. 
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FIG. 1: Phase diagram predicted from MF calculations for the cell model with fixed H bond 
strength (J = 0.5), fixed H bond volume increase (vhb/'I'o = 0-5), and different values of the H 
bond cooperativity strength J„. (a) Singularity-free scenario (Jo- = 0). At high T, liquid (L) 
and gas (G) phases are separated by a first order transition line (thick line) ending at a critical 
point C, from which a L-G Widom line (double-dot-dashed line) emanates. In the liquid phase, 
maxima of Kt and ap converge to a locus (dot-dashed line). At C both Kt and ap and have 
diverging maxima. The locus of the maxima is related to the L-L Widom line for Tc — > (see 
text), (b) Liquid-liquid critical point scenario with positive critical pressure (for Ju = 0.05). At 
low T and high P, a high density liquid (HDL) and a low density liquid (LDL) are separated 
by a first order transition line (thick line with HDL/LDL labeled) ending in a critical point C", 
from which the L-L Widom line (dot-dashed line) emanates. Other symbols are as in the previous 
panel, (c) Liquid-liquid critical point scenario with negative critical pressure (for J„ = 0.35). Here 
the L-L Widom line (dot-dashed line) is shown intersecting the L-G spinodal (dotted line). Other 
symbols are as in the previous panel, (d) Critical-point free scenario (Jo- = 0.5). The HDL- 
LDL coexistence line extends to the superheated liquid region at P < 0, reaching with the liquid 
spinodal (dotted line). The stability limit (SL) of water at ambient conditions (HDL) is delimited 
by the superheated liquid-to-gas spinodal and the supercooled HDL-to-LDL spinodal (dashed 
line), giving a re-entrant behavior as hypothesized in the SL scenario. Other symbols are as in the 
previous panels. In all panels, is the Boltzmann constant. 
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FIG. 2: Possible scenarios for water for different values of H bond energies J, (directional com- 
ponent) and Jo- (cooperative component), both in units of the van der Waals energy e, obtained 
from MF calculations. The ratio vhb/vq is kept constant, (i) If Jo- = (red line along x-axis), the 
singularity free (SF) scenario is realized, independent of J. (ii) For large enough J^-, water would 
possess a first-order liquid-liquid phase transition line terminating at the liquid-gas spinodal — the 
critical point free (CPF) scenario; the liquid spinodal would retrace at negative pressure, as in the 
stability limit (SL) scenario (yellow region in top left), (iii) For other combinations of J and J^, 
water would be described by the liquid-liquid critical point (LLCP) scenario. For larger J^-, the 
LLCP is at negative pressure (ochre region between dashed lines). For smaller Jo, the LLCP is 
at positive pressure (orange region in bottom right). Dashed lines separating the three different 
regions correspond to mean field results of the microscopic cell model. Equations for the lines 
are Jo = J/2 and Ju = a + bJ, with a ~ 0.3 and b ~ 0.36. The P — T phase diagram evolves 
continuously as J and Jo change. 
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FIG. 3: Phase diagrams from MC simulations. We fix the parameters J = 0.5, uhb/^o = 0.5 and 
vary the parameter J^-. (a) For J^- = 0.5 (i.e. > a + hJ where a and h are given in the text and 
both are approximately 1/3), we find a liquid-liquid phase transition (thick line with circles) ending 
on the superheated liquid-to-gas spinodal (dotted line with diamonds) as in the CPF scenario, (b) 
For Jo- = 0.3 (i.e. > J/2), the liquid-liquid phase transition ends in a liquid-liquid critical 
point (LLCP) at negative pressure, (c) For Jg- = 0.05 (i.e. J^- < J), the LLCP ends at positive 
pressure and the line of specific heat maxima (crosses), emanating from the LLCP, is shown only 
for positive pressure. Errors are of the order of the symbol sizes. Lines are guides for the eyes. 
Other model parameters are as for MF calculations (see text). 
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